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1.  INTRODUCTION 

It  has  long  been  appreciated  that  the  simple  blrth-and-death  Markov 
process  often  provides  an  adequate  Initial  model  for  the  behavior  of 
service  systems,  populations,  epidemics,  and  many  other  stochastic 
systems;  see  Feller  [3]  for  an  early  classic  account.  Refinements, 
particularly  in  the  modeling  of  service  systems,  have  typically  Involved 
the  replacement  of  exponential  service  times,  also  assumed  Independent,  by 
independent  random  variables  of  "general"  or  arbitrary  distribution, 
replacement  of  independent  exponential  inter-arrival  times  by  independent 
"general"  random  variables,  or  both.  Isolated  examples  in  which  arrival 
process  parameters  are  allowed  to  change  deterministically  in  time  have 
also  been  studied. 

Randomly  appearing  fluctuations  in  system  environment,  associated 
with  weather  or  other  change  in  physical  surroundings,  personnel  changes, 
alteration  of  system  usage  intensity,  etc.,  are  likely  to  be  reflected 
in  changes  in  observed  failure  and  repair  rates  ir.  a  system  reliability 
context. 

1.1  Examples  :  In  order  to  illustrate  the  ideas  summarized  above,  we 
introduce  two  specific  birth  and  death  models  in  random  environments. 
Numerical  results  and  discussion  will  be  found  in  the  last  section. 

Repairman  Model.  Suppose  any  of  m  machines  that  are  in  use  fail  independent¬ 
ly  in  Markovian  fashion  at  rate  Ut>  t  >  0),  and,  if  failed,  experience 
repair  at  rate  (ut,  t  >  0}.  In  turn,  the  rates  Ut)  and  (ut>  are  them¬ 
selves  finite-state  Markov  processes  independent  of  the  state  of  the 
machines,  so  that  if  n  is  the  number  of  machines  on  repair,  and  j  identifies 
the  environmental  state,  then  (n,j)  is  the  overall  state  variable  of  the 


system;  conditional  upon  j,  n  cnanges  Dy  one  unit  at  a  time  in  typical 
birth-death  fashion  with  parameters  depending  on  j  :  the  probability  that 
an  operating  machine  goes  down  in  (t,t+dt)  is  x^m-njdt  +  o(dt),  while  tne 
probability  that  a  machine  on  repair  becomes  available  is  Wj  min(R,n)dt  +  o(dt)  , 
R  being  the  number  of  repairmen.  Such  a  setup  describes  groups  of 
redundant  equipments  that  all  experience  common  environmental  intensities 
simultaneously;  the  environmental  changes  are  reflected  in  the  numerical 
values  of  the  failure  and  repair  rates  that  prevail  at  any  time  point. 

The  model  described  also  may  be  used  to  represent  the  behavior  of  a  net¬ 
work  of  m  timesharing  computer  terminals  that  independently  send  messages  or 
programs  to  a  common  central  computing  facility.  Suppose  the  rate  of  message 
transmission  changes  with  external  activity,  i.e.  responds  to  an  occasional 
period  of  unusual  activity,  perhaps  a  crisis  situation.  In  tnis  case  the 
"constant"  demand  rate  X  (assumed  equal  for  all  terminals)  switches  almost 
instantaneously  to  a  higher  value,  switching  back  to  normal  after  the  crisis 
elapses,  but  doing  so  repeatedly.  Left  alone,  the  central  processor  may  well 
continue  processing  at  the  original  rate,  allowing  congestion  to  simply  build 
up,  and  eventually  drop  again  when  the  demand  lapses.  On  the  other  hand, 
remedial  action  may  be  taken.  These  phenomena  suggest  interesting  and  realis¬ 
tic  questions  concerning  stochastic  control,  but  these  will  not  be 
considered  in  this  paper. 

Mass  Search  Model .  Let  a  group  of  m  predators  attempt  to  round  up  and  capture 
a  finite  group  of  p  prey.  Predators  search  independently  and  with  equal  inten¬ 
sity  and  effectiveness,  but  the  detection  rate  is  allowed  to  depend  upon  exter¬ 
nal  environmental  conditions  that  cause  relatively  long-term  changes  in,  say, 
visibility.  If  a  prey  is  detected,  it  is  followed  until  lost  by  the  predator, 
after  which  moment  it  is  susceptible  again  to  search,  detection,  and  active 
surveillance.  While  predators  and  prey  move  independently,  as  soon  as  a  pre¬ 
dator  begins  following  a  particular  prey  the  latter  is  removed  from  circula¬ 
tion  and  the  remaining  unattached  predators  continue  search  for  the  free  prey. 
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Of  interest  is  the  long-run  or  stationary  distribution  of  the  number  of 
prey  under  simultaneous  surveillance,  and  also  the  first-passage  time  until 
a  large  fraction-  perhaps  all-  of  the  prey  are  simultaneously  under  the  eye 
of  the  predators.  Here  is  a  plausible  model  :  the  state  of  the  system  is 
(n,j),  so  that  if  the  number  of  prey  under  surveillance,  is  n,  and  the 
state  of  the  environment  is  j,  then  the  probability  that  a  free  prey  is 
detected  is  Xj(p-n)(m-n)dt  +  o(dt),  while  the  probability  that  a  prey  under 
surveillance  is  lost  is  pjndt  +  o(dt),  where  n  <  m1n(p,m)  and  Xj  and  Uj 
change  in  accordance  with  Independent  finite-state  Markov  processes,  as 
before . 

1.2  The  analytic  structure  :  Me  shall  consider  Markov  processes 
{(Xt,Yt),  t  >  0),  on  the  state  space  {(n,j),  0  <  n  <  N,  1  <  j  <  Kn> ,  with  a 
block-tridiagonal  infinitesimal  generator  Q  : 
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Here  4 


a^e  square  matrices,  respectively  of  order 


An- . a(n) 

,  Kj,  .  Their  diagonal  elements  are  strictly  negative,  all  other 

elements  are  non-negative.  The  matrices  A^,  0  <  n<  N-l,  and 
1  <N,  are  rectangular,  with  appropriate  dimensions  ;  their  entries 

are  non-negative.  The  rowsums  of  Q  are  equal  to  zero,  therefore  we 

have  that 

A;o)e  +  A(o)e  x  0  > 

M^)e  +  A^e  +  =  0,  1  <  i  <  N-l, 

M^e  +■  A^e  =  0, 

wher*3  e  denotes  a  column  vector  with  unit  elements.  The  variable 
is  to  be  interpreted  as  the  state  of  the  environment,  and  X  as  the 
state  of  the  birth  and  death  process,  at  time  t  . 

We  assume,  furtnernore,  that  the  Markov  process  Q  is  irreducible,  and  we 
denote  by  level  n  the  set  {(n,j),  1  <  j  <  Kn),  of  states  corresponding  to 
the  common  value  r.  for  the  first  index.  The  structure  (1.1)  of  Q  permits 
the  Markov  process  to  move  up  or  down  by  only  one  level  at  a  time.  It  is 
in  this  respect  analogous  to  the  classical  birth-and-death  process,  see 
Feller  [3]  or  Karlin  and  Taylor  [i]. 

In  the  examples  cited  above,  either  X  or  Y  changes  each  time  the 
Markov  process  undergoes  a  change  of  state.  In  fact,  this  restriction  is 
rot  part  of  the  model,  and  we  allow  X  and  Y  to  change  simultaneously. 

1.3  txi  sting  literature:  Infinite  birth-and-death  models  in  random 
environment  have  been  studied  for  some  time  already;  early  results  are 
found  in  Eisen  nd  Tainiter  [2],  Purdue  [11]  and  Yechiali  [14].  More 
recently,  Neats  [9],  Chapter  6,  has  systematically  examined  a  class  of 
problems  in  which  N  =  and  Q  has  a  special  repetitive  structure  : 


This  structure  leads  to  matrix-geometric  stationary  probability  vectors  and 
efficient  algorithmic  procedures. 


5. 


Finite  models  have  been  examined  by  Torrez  [12],  who  suggests  using 
numerical  procedures  designed  to  solve  eigenvectors  for  band-matrices. 

Hajek  has  considered  in  [4],  Section  5,  a  finite  model  with  repetitive 
structure  : 

-A<2>  -  ...  -  A<N-2>, 

Ad)  =a(2)  .  ...  .A(H-l)t 

m{2)  =m(3)  „  ...  .M(M-l)f 

K1  =  *2  =  •••  = 

Hajek  determines  the  stationary  probability  vector  in  terms  of  two  matrices 
R  and  R,  of  order  K,  which  have  to  be  iteratively  computed.  Finally,  Keilson 
et  al.  [6]  have  considered  finite  models  with  the  structure  (1.1)  and  equal 
Kn's.  They  analyse  the  (Laplace  transform  of)  first  passage  times  distributions, 
from  which  they  obtain  equations  for  moments,  and  for  the  stationary  distribu¬ 
tion.  We  cannot  in  this  short  space  describe  in  detail  the  differences  between 
our  results  and  those  in  [4,  6,  12];  but  we  shall  give  some  discussion  at 
the  end  of  Section  5. 

This  paper  presents  an  efficient  computational  approach  to  the  analysis 
of  birth-and-death  models  in  a  Markovian  environment.  The  emphasis  is  upon 
obtaining  numerical  properties  of  both  stationary  distributions  (in  the  next 
section),  and  first-passage  time  (;n  Sections  3  and  4).  The  computational 
algorithms  are  discussed,  and  numerical  examples  are  given,  in  the  last  two 
sections. 

2.  THE  STATIONARY  DISTRIBUTION 

In  order  to  determine  the  stationary  probability  distribution,  and  moments  of 
first  passage  times,  it  is  useful  to  think  of  the  Markov  process  as  evolving  in 
a  certain  manner. 


6. 


Fur  0  •-  n  v  define  Sn  to  be  the  restriction  of  the  original  process  Q, 
oDserveo  during  tnose  intervals  of  time  spent  at  level  n,  before  the  original 
process  enters  level  n+1  for  the  first  time,  the  state  space  of  Sn  is 
<(n,j),  1  j  Kni.  Clearly  all  S  ,  0  <  n  <  N-l  are  transient  Markov  proces¬ 
ses.  Toe  process  S,,  is  the  restriction  of  the  process  Q  to  the  states 
{(N ,j } ,  1  h  j  s  K  };  it  is  an  ergodic  Markov  process.  We  denote  by  the 
infinitesimal  generator  of  the  process  Sn,  0  <  n  <  N. 

In  order  to  determine  the  matrices  C  ,  we  need  the  following  result. 

Le—i  1 

Consider  a  varkov  process  on  the  state  space  (1,2, . . . ,r,r+l,r+2 . r+sl, 

w.n  infinitesimal  generator  Q  : 


w^ere  A  is  a  square  matrix  of  order  r,  A  is  a  rectangular  r  by  s 
"■  i •  r i * .  T>e  states  r+1  to  r+s  are  all  absorbing. 

>  nhp  ^sates  1  to  r  are  all  transient  if  and  only  if  the  matrix  A  is  non- 

Sinqu’ar. 

.  'Ke  { i , i ) th  entry  of  {-A'1),  for  1  <  i.j  <  r;  is  the  expected  amount  of 
time  soeot  in  the  transient  state  j,  starting  from  the  transient  state  i, 
before  jbsorotion  in  any  of  the  absorbing  states. 

,  Toe  (i,k)t.h  entry  of  (-A  *  a),  for  1  <  i  <  r,  r+1  <  k<r+s,  is  the 
probability  that,  starting  from  the  transient  state  i,  absorption  occurs 


in  the  state  k. 
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Proof.  The  first  assertion  is  proved  in  deuts  [9],  Lemma  2.2.1. 

The  matrix  A**  has  all  nonpositive  entries.  The  second  assertion  is 
a  consequence  of  Neuts  and  Meier  [10],  Corollary  2.  To  determine  the 
probability  P^  k  of  being  eventually  absorbed  in  state  k,  r+1  <  k  <  r+s, 
starting  from  state  i,  1  <  i  <  r,  we  study  the  Markov  chain  embedded  at 
instants  immediately  following  a  transition  in  the  Markov  process. 

'V 

The  corresponding  transition  probability  matrix  P  is  given  by 


p  = 


I  -A  aA 


-a"!A 


0  I 


where  the  matrix  A  is  diagonal,  with  diagonal  entries  equal  to  those  of  A. 
It  results  from  Kemeny  and  Snell  [7],  Theorem  3.3.7,  p.  52,  that 

P  =  [I  -  (I  -A^A)]"1  ( -A_1A)  =  ( -A- 1 ) A ; 
which  completes  the  proof  of  Lenma  1.  □ 


We  may  now  determine  the  matrices  Cn. 
Lemma  2 


C  =  A(")  +  M<n)(-C:11)A(n‘1).  1  <n  <N. 
n  n  '  * 

Proof.  The  equation  (2.2)  is  obvious  starting  from  level  0,  the  process  SQ 
terminates  as  soon  as  the  process  Q  enters  level  1.  Meanwhile,  the  transi¬ 
tions  are  governed  by  A(0).  Since  the  process  Q  is  irreducible,  the 
process  SQ  contains  only  transient  states,  and  the  matrix  CQ  is  non¬ 
singular,  by  Lemma  1.  Also,  (-Cg)  >0. 


(2.2) 

(2.3) 


the  generator  of  the  process  S^.  Let  Z^(t), 
fcilews.  Z,(x)  =  i  if  the  process  S1  is  in  state 
site  that  Z,(t)  =  i,  and  Z,(x  +  dr)  =  j  *  i. 

x  1 


;  i  tieh  l  ast  have  occurred  in  the  process  Q.  Z^  (a  +  d-t)  =  j 
f  the  following  events  have  occured. 
f-r.i  ;i,i,  to  (l,j),  this  happens  with  probability 


- !  ; .  to  (0,k),  for  some  k,  and  the  process  Q 

after  spending  an  unspecified  amount  of  time  at  level  0, 
y  ot^er  state  at  level  1.  This  happens  with  probability 
(dx),  where  y.  .  is  the  probability  of  moving  from 
ore  visiting  any  other  state  at  level  1. 

a  =  A-0'  and  A  =  it  results  from  Lemma  1  that 


I'Zjft)  »  i]  «  Aj1]  dx  +  E°  Mj1)  (-A(0)_1A(0))k  i 


dx 


+  o  (dx). 


1  <  i*  j  <  Kj  , 


1  Z !  ( t 1  =  i] 


(1+aJ1)  dx)  + 

'  t  1 


K<‘>  (-aW-V")),  ,  dT 


k=l 


i  ,k 


'k,i 


1  • 


+  o  (dx). 


1  <  i  <  K 


9. 


Thus,  Cj  =  A^U  M^f-A^rV0), 

.  a<1)+ 

which  proves  (2.3)  for  n*l. 


Assume  that  (2.3)  holds  for  n,  we  prove  now  that  it  holds  for  n+1.  Let 
Zn+i  (t),  t  >  0,  be  equal  to  i  if  the  process  Sn+1  Is  In  state  (n+l,i) 
at  time  t. 

If  Zn+1  (t)  =  i,  then  Zn+1  (r+dr)  *  j*i  if  and  only  if  one  of  the  fol¬ 
lowing  events  occurs. 

a.  There  is  a  transition  from  (n+l,i)  to  (n+l,j),  with  probability 
A^t ^dr  +  o(dx) . 

b.  There  is  a  transition  from  (n+1,1)  to  (n,k)  for  some  k,  with  probability 
M^j^dt  +  o(dx),  and  the  process  Q  returns  to  (n+l,j),  after  spending  an 
unspecified  amount  of  time  at  levels  0,1,..., n,  before  visiting  any  other 
states  at  level  n+1,  with  probability 


Since  Sn  records  the  visits  made  by  Q  at  level  n  before  reaching  level 
n+1,  we  have  by  lenma  1  that 


,(n)i 


(2.4) 


It  is  now  easy  to  prove  that  (2.3)  holds  for  n+1,  which  completes  the  proof 
of  the  lemma. 


o 


i"  the  stationary  probability  vector  P,  i  .e . 
r  ^  i  =  We  Partition  that  vector  as 
v  the  subvectors  have  Kn  elements  and 
.  ::f  I*- vel  n,  0  <  n  <  N. 


’’  '  N,  are  determined  by  the  equations 


.-1 


) ,  for  0  <  n  <  N-l, 


(2.5) 

(2-6) 


(2.7) 


.l*  j,-<?  CM)  of  Q,  it  results  that  the  system 

•■■''ose d  into 


c-i 


+  P  M^n+1^ 
-n+1 


M  -  0. 


(2.8) 

=  0,  1  <  n  <  N-l,  (2.9) 

(2.»! 


■->  are  matrix  equivalents  of  the  familiar 

for  birth-and-death  processes.  They  may  also  be 


the 


;;r'»l',ihility  balance. 
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From  (2.9),  Pg  =  -Pj  ^  )'1=  Pj  M^^-Cg1).  Then,  by  recurrence, 

using  Lemma  2,  the  vectors  P^  ,  1  <  n  <  N,  satisfy  the  equations  (2.5)  and 
(2.6).  Since  the  matrix  is  the  infinitesimal  generator  of  a  finite 
irreducible  Markov  process,  Equation  (2.5)  has  a  unique  solution,  up  to  a 
multiplicative  constant,  and  that  constant  is  determined  by  (2.7). 

□ 

This  result  suggests  an  algorithm  to  compute  the  stationary  probability 
vector. 

Algorithm  A. 

Al.  Determine  recursively  the  matrices  Cn,  0  <  n  <  N. 

A2.  Solve  the  system  =  0,  e  =  1. 

A3.  Compute  recursively  the  vectors  P^  ,  n  =  N-l,  ...,  0,  using 
instead  of  P^  . 

A4.  Re-normalize  the  vector  P  so  obtained. 

A  complete  analysis  of  this  algorithm  is  deferred.  We  shall  make  some 
corrments  on  it  at  the  end  of  this  section.  Examples  of  numerical  appli¬ 
cations  to  certain  specific  models  are  presented  at  the  end  of  the  paper. 

3.  FIRST  PASSAGE  TIME  TO  HIGHER  LEVELS 

We  denote  by  Tn  the  first  passage  time  from  level  n-l  to  level  n, 
and  by  xn  m  the  first  passage  time  from  level  n-l  to  level  m  >  n  : 

m  =  inf  { t  >  0  :  =  m |  Xg  =  n-l)  ,  1  <  n  <  m  <  N, 


We  define,  for  x  >  0,  1  <  n  <  N,  1  <  i  <  j,  1  <  j  <  Kn> 


nin!(x)  =  P  ET  <  x,  Yt 


=  j  I  X„  =  n-l,  Y„  =  il. 


c  state  of  the  system  at  time  h,  a  simple  probabilistic 
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a+A^n”1)h)g<nhx)  +  I  A^j^h  gj;nj(x)  + 

1,1  1,J  k*j  l.K  K,J  1,J 


*  ifr"'  ife"'  •£!“>  * 0(h)’ 


(3.1) 


denotes  the  Stieltjes  convolution.  Substracting 
•les  of  equation  (3.1),  dividing  by  h  and  letting 
first  order  differential  equation  for  g(nhx). 

•  *J 

,  toe  Laplace-Stieltijes  transform  of  gjn|(x), 


g(n}(x)]  dx,  and  by  Gl";(£)  the  matrix  with  entries 

0  X  1  „  J 


;(n) 


j.l)  that 

i  =A(n'1Mn)(OtA(n-l)  +  f|<n’1Mn"1J(5)G<nJ(c). 

(3.2) 

for  2  <  n  <  N, 


-  Al°>S(1)(5)  .  A<°>.  (3.3) 

ho  written  as 

90(O  A(0),  (3.4) 

-  D  j(0  A(n’1),  for2<n<N,  (3.5) 

(M  -  (3.6) 

-  A(n)-  M(nV,(n)(0)-1,  1  <  n  <  N-l .  (3.7) 
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Ut  -Mr  .<«.  »,  ■  j  I  <0  .  n-1.  V0  .  )]  . 

n  ,m 

and  let  G(n;m^(5)  denote  the  Laplace-Stieltjes  transform  of  g(nlm^(x). 

'  * J  1 tj 

We  readily  obtain  that 

G(n,m)(5)  .  G(n,m-1)^  G^m^(5), 

=  M  G(k)(0. 

k=n 

for  1  <  n  <  m  <  N,  where  we  define  G^n,n"^(5)  =  I  for  all  n. 


(3.8) 

(3.9) 


We  easily  prove  from  Lemma  2  and  Equations  (3.4)  to  (3.7)  that 

°n  =  -  Cn1  »  for  0  <  n  <  N-1,  (3.10) 

hence  G(n)(0)  =  •  A1"'1*,  for  1  <  n  <  N,  (3.11) 

which  is  merely  another  representation  of  Equation  (2.4}  ,  since  (3^(0)).  . 
is  the  probability  that,  starting  from  (n-1,1),  the  process  0  visits 
(n,j)  before  visiting  any  other  state  at  level  n. 

Let  U(n’m>  =  -  —  G^’^U)],  n  . 

35  ^ 


We  have  that 


Jj"f>  =  E[t  Y 
i  ,J  n  ,m  t 


+0 


=  j  |  Xn  =  n-1,  =  i ]  , 


n,m 


*"d  -  E<'n,„l*0  .  „-l,  v0  •  I]  . 

where  u^n,m^  =  U^n,m^e. 

We  define  =  u^n,n^  - - G^(5)|  n  ,  and  u^  =  U^e 

35  5 


14. 


:.s  of  the  first  passage  time  to  higher  levels  satisfy 
■:e  relations. 


r-1  (e+M<n‘1)u(n'1)),  for  2  <  n  <  N, 

“n-1  -  - 

+  G(n,m'1)(0)  u(m\  for  1  <  n  <  m  <  N. 


(3.12) 

(3.13) 

(3.14) 


\j  be  proved  by  differentiating  equations  (3.6)  and 

and  evaluating  the  derivatives  gjr  D  (f.)  at  i  =  0  . 
■ >  on  first  probabilistic  principles.  First,  equation 
t.-mma  1  (b)  and  (2.2).  Then  let  V  be  the  elapsed 
■,  enters  a  different  level;  that  is, 

X,j '  .  By  the  strong  Markov  property,  Lemma  1,  and 


Y0 


X0  n  ’  Y0  1  ^ 


n  -  1  [X  =  n  ,  Yn  =  i] 


n  - 1 


k  - 1 


A(n))-1  M(n)Jik  E[TJX0  =  n  -  1  ,  Yft  =  k] 


M(n}  iLt-C-^)  A^n-^:)ijE[Tn+1|X0  =  n  ,  YQ  =  j] 


.15) 


15. 


Hence, 

[I  .  M(n)(-C^,)A('>-,>]  u(r,+1) 

I 

=  (-a^V1  e  +  (-A(n))"1  M(n)  u(n) 

Multiplying  both  sides  of  (3.16)  by  results  in  the  equation 

[A(n)  +  u(n+1) 

=  -[e  +  M<n>  u(n)J. 

Equation  (3.13)  now  follows  from  (2.3)  and  (3.17). 

Since 

E£VmlX0'"-'-Vi:] 

=  u<"'n,-,)(i)  *  G(n'"-,)(0)  M  (,■) 
equation  (3.14)  also  follows. 


(3.16) 


(3.17) 


Any  number  of  moments  of  the  first  passage  times  may  be  similarly  obtained. 
We  merely  state  the  following  result  for  second  moments,  without  proof. 


Example:  Consider  the  repairman  model  of  Section  (1.1).  Assume  there 
are  just  two  environment  states,  denote  by  j  =  1,  2  .  Let  the  transition 
rate  from  environment  state  2  to  1  be  a  ,  and  from  environment  state 
1  to  2  be  £  .  Then  G^(c)  satisfy  the  following  system  of  equations. 


A,(m) 


.0)  /,x  _  Ar'"'  „  /1X  .  a  „(1) 

Gl,j  {c)  "  d^l)  +  5  V1}  d-|  (1 )  +  1  G2,j  U) 


4!]  (c)  =  d2(0  +  l  £j(2)  +  d2 ( 1  f  +  I  G1J  U) 


A9(m) 


6 _ ,0) 


(3.21) 


c(n+l)  -  Xl(m‘n)  t  m  ,  ulCmin(R*n)]  r  (n)  ,  x  c(n+l)  ,  n 
Gl,j  U)  -  dj(n TV l  *j(1)  +  d^n)  LG  U)  G  U)]lp- 


.  a  ^(n+1) 

*  d, (n)  +  5  G2,j  (t)  • 


(3.22) 


c(n+l)  M  _  X2(rn‘n)  {0]  u2[min(R,n)]  (n)  (n+l) 

G2,j  “  d2THl’TT  *j(2)  +  d^n)  +1  ~  [G  G  U)]2,j 


,  6  r(n+l)  /  , 

+  S^)  —  G1  ,j 


for  n  <_  m  -  1 
where 


d^l) 

=  X-jm  +  a 

» 

d2(l) 

=  A2m  +  8 

) 

d-|  (n) 

=  X-j(m-n) 

+  y-|[min(R,n)]  +  a  ; 

d2(n) 

=  x2(m-n) 

+  u2[min(R,n)]  +  6  ; 

fl 

if  i  =  j  , 

M1) 

- 

(3.23) 


otherwise 
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Vie  above  equations  can  in  principal  be  solved  recursively.  Finally, 

G0  »n+1 )  =  g(1)(0  G^U)  x...x  G(n+1)U)  • 

Similarly,  the  expected  first  passage  times  satisfy  the  following 
recursive  equations. 


(3.24) 


o  =  -  dTOT  +  H-frT 


T  I  Y 
Jl'r 


E  Ti  iYo " 2]  ■  d^nr +  d^fry  ECTi  lYo = 1] 


nrliY0  ^ 


dyrnT  +  d^THT  ECTn+l  (Y0  =  2] 

u,[min(R,n)]  f  -» 

+  -Ld70r5 - HlYo  =  1^WYo  =  1l[ 


V,iyo-2’sdjUtd7HTE[Vllvo  =  ,] 


M9[min(R,n)]  r  -> 

1  lcrT  =  2]  +  E[T„X1|V0  =  2]] 


+  d9(n) 


{«T„I»0  = 


n+1 
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4.  FIRST  PASSAGE  TIME  TO  LOWER  LEVELS 

There  is  some  symmetry  in  the  process  Q,  which  we  have  not  exploited 
yet.  This  we  proceed  to  do  now.  Instead  of  the  processes  Sn,  0  <  n  <  N, 
defined  earlier,  we  now  consider  the  processes  Sn,  0  <  n  <  N. 

For  1  <  n  <  N,  Sn  is  the  restriction  of  the  process  Q,  observed  during 
those  intervals  of  time  spent  at  level  n,  before  the  process  Q  moves  down 
to  level  n-1  for  the  first  time.  All  §n,  1  <  n  <  N,  are  transient  Markov 
processes.  The  process  5g  is  the  restriction  of  Q  observed  at  the  lowest 
level;  it  is  an  ergodic  Markov  process.  We  denote  by  Cn  the  infinitesimal 
generator  of  the  process  Sn>  0  <  n  <  N. 

We  may  carry  out  for  the  processes 5n  exactly  the  same  analysis  as  we 
did  for  the  processes  Sn.  We  indicate  below  the  main  results  for  two 
reasons . 

a.  The  first  passage  times  to  lower  levels  must  be  analysed  via  the 
processes  $n> 

b.  We  shall  be  able  to  describe  the  precise  correspondence  between  our 
analysis,  and  the  analysis  of  Neuts  for  the  infinite  quasi-birth-and- 
death  process. 

The  proof  of  the  next  lemma  is  omitted,  since  it  is  identical  to  that 


of  Lemma  2  and  Theorem  1. 


iO  . 


Lemma  3 . 

The  matrices  Cn  ,  0  <  n  <  N,  are  recursively  determined  as  follows. 


=  A(n)+  A(n)(-Cn+1_1)  M(n+1)  ,  for  0  <  n  <  N-l. 

The  stationary  probability  vector  P_  =  (Pg,  P^,  ....  P^)  is  determined  by 
the  equations 


Po£o  =  -  ’ 


^  ■  Vl  A(n_1,(-C-1),  for  1  <  n  <  N, 


N 

I 

n=0 


Pn 


e  =  1. 


(4.1) 

(4.2) 


□ 


The  (i,j)th  entry  of  the  matrix  A^n  *^(-C  *)  is  equal  to 


[A 


K  A(n-1) 

*n  Ai,k 


(4.3) 


-  MjV))  1  —  (<)kj  • 

1,1  k=l  a(n-l)  n  K,J 
i 

The  factor  A^nJ1V(-A|nT^)  is  the  probability  that,  upon  leaving  the  state 
(n-l,i),  the  process  Q  moves  to  the  state  (n,k).  The  (k,j)th  entry  of 
(-C'1)  is  the  expected  time  spent  by  the  process  in  state  (n,j),  starting 
from  (n,k),  before  hitting  any  state  at  level  n-l  (by  Lemma  1).  Therefore, 
it  results  from  (4.3)  that  [A^n_1^(-C"1)]i  j  is  equal  to  (-a|”T^)  times 
the  expected  time  spent  in  the  state  (n,j),  before  the  first  return  to  level 
n-l,  given  that  the  process  9  starts  in  the  state  (n-l,i).  This  is  exactly 
the  interpretation  of  Neuts'  matrix  R  for  the  infinite  quasi-birth  -and- 


1*. 


Let  t  „  denote  the  first  passage  time  to  level  n,  down  from  level 
m,n  3 

m+1,  for  0  <  n  <  m  <  N-l  : 


Vn  =  inf  {t>0  :  Xt  =  n  I  X0  =  m+1}> 


and  let 


E  ‘Vn  I  *0  '  «»•  V0 

d^  =  d^’0}  , 

1  *»  * 


=  11  , 


for  0<n<m<N-l.  The  proof  of  the  next  theorem  is  omitted,  since  it 
is  identical  to  that  of  Theorem  2. 


Theorem  3. 

The  expected  values  of  the  first  passage  time  to  lower  levels  satisfy  the 
following  recurrence  relations. 

d(N_1)=  -  e  ,  (4.4) 

d(m)=  (j(m»n)=  -  C“jj(e  +  A(m+1)d(m+1>),  for  0  <  m  <  N-2,  (4.5) 

d(m’n)=  d(m’n+1>  +  g(m»n+1) (0)^(n)  ,  for  o  <  n  <  m  <  N-l.  (4.6) 

where  S(m,n)(0)  =  n  (-Ck+1)  M(k+1> ,  for  0  <  n  <  m  <  N-l,  (4.7) 

k=m 


Observe  that  in  the  right-hand  side  of  (4.7)  »  the  left  most  matrix 
in  the  product  is  not  (-  C^). 
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5.  REMARKS  ABOUT  THE  COMPUTATIONAL  ALGORITHMS 

The  algorithm  A  described  in  Section  2  is  numerically  stable 
under  a  large  range  of  values  for  the  entries  of  the  generator  Q 
of  (1.1).  The  matrices  (-C  *)  have  only  non-negative  entries.  Moreover, 

since  they  measure  the  time  spent  at  level  n  only,  before  moving  to 
level  n+1,  they  usually  are  of  the  same  order  of  magnitude  for  each  n, 

-;,'en  for  large  values  of  N,  if  the  entries  of  the  matrices  are  of 
toe  same  order  of  magnitude  for  each  n,  and  similarly  for  the  matrices 
a  potential  source  of  trouble  exists  when  the  A's  are  either  very  small  or 
v 9 '"j  large,  compared  to  the  other  elements  of  Q.  In  that  case,  the  expected 
~  -nes  spent  at  one  level  before  moving  up  are  respectively  very  large  or 
very  small,  and  there  is  a  risk  of  encountering  overflow  or  underflow 
problems,  when  determining  the  matrices  (-C”*). 

The  steps  A3  and  A4  of  Algorithm  A  are  more  delicate.  We  start  step  A3 
vH-h  a  vector  tt^  normalized  by  e  =  1.  If  N  is  very  large,  it  is  likely 
to at  the  vector  will  be  much  smaller  than  ,  and  there  is  a  real  risk 
of  running  into  overflow  problems  while  performing  step  A3.  In  order  to 
;  vercome  this  difficulty,  we  have  merged  the  two  steps  A3  and  A4,  and 
•■-normalized  the  vectors  each  time  a  new  subvector  is  determined  (see 
Algorithm  B  in  Appendix  A). 

It  results  from  Theorem  2  that  one  may  compute  the  vectors  u^n,m^  at 
the  same  time  as  one  is  preparing  the  evaluation  of  the  stationary  probabi¬ 
lity  distributions,  i.e.  during  step  A1  of  Algorithm  A.  In  the  numerical 
pv’mples  presented  in  the  next  section,  we  have  determined  the  first  passage 
+  imos  from  level  0  to  level  m,  for  1  <m  <  N.  The  corresponding  algorithm 


.  a 
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Finally,  we  compare  the  numerical  efficiency  of  three  approaches. 

As  observed  by  Torrez  [12],  the  matrix  Q  is  a  band-matrix,  and  there  exist 

efficient  numerical  procedures  to  solve  eigenvectors  for  band-matrices. 

'S 

The  complexity  of  such  procedures,  for  the  matrix  Q,  is  0(N,K  ),  where 

K  =  max{Kn,  0  <  n  <  N)  (Wilkinson  and  Reinsh  [13],  p.  70).  The  crucial 

step  in  our  algorithm  resides  in  the  inversion  of  the  matrices  Cn,  0  <  n  <  N-l, 

and  the  solution  of  the  system  (2.5). 

M  3 

The  corresponding  complexity  is  0(  I  K  ).  In  Keilson  et  al .  [6],  2  N  matrices 

n=0  n 

have  to  be  inverted,  and  N  systems  have  to  be  solved.  Again,  the  complexity 
N  3 

is  0(  I  K  ) .  Clearly,  the  three  approaches  have  globally  similar  numerical 
n=0  n 

efficiency;  in  order  to  distinguish  among  them,  one  would  have  to  determine  the 
coefficients  implicit  in  the  0(.)  notation.  However,  it  must  be  observed  that 
the  last  two  have  the  additional  advantage  of  offering  clear  probabilistic 
interpretations  of  the  computed  quantities. 

6.  NUMERICAL  RESULTS 

6.1  Machine  repairman  model.  We  consider  a  system  with  N=5  machines  and 
R=1  repairman.  The  system  can  be  in  K=2  environments.  The  system  remains 
in  the  environment  state  j  (j  ■  1,2)  for  an  exponentially  distributed  random 
interval  of  time,  with  parameter  c^.  The  failure  rate  of  each  machine  is 
equal  to  Xj  in  the  environment,  with  x^  =  0.12,  and  X£  *  0.06.  The  repair 
rate  of  the  repairman  is  u  =  1  in  both  environments. 

Our  objective  is  to  measure  how  the  rate  of  changes  in  the  environment 
influences  the  system  behaviour.  In  order  to  do  so,  we  set  *  a/2,  ■  B, 

and  chose  different  values  for  a.  Then  the  stationary  probabilities  of  being 
in  each  environment  remain  constant  and  are  given  by  yj  =  2/3  and  >2  a  1/3. 

We  expect  that  if  e  is  large,  then  the  environment  changes  rapidly,  and  the 
system  is  only  influenced  by  the  average  failure  rate  XQ  *  >2x2  s  °*1* 
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On  the  other  hand,  if  g  is  small,  then  the  environment  stays  for  long 
periods  of  time  in  the  same  state,  and  this  should  affect  the  dynamic 
de-navi  our  of  the  system. 

We  denote  by  c^(e),  0  <  i  <  5,  the  marginal  distribution  of  the 

number  of  machines  on  repair,  for  a  given  value  of  6. 

We  furthermore  denote  by  n^(X),  0  <  i  <  5,  the  probability  distribution 

for  the  classical  machine  repairman  system,  with  constant  failure  rate  X. 

In  table  I,  we  give  the  cumulative  probabilities,  corresponding  to 

-5  -2  -1 

for  x  =  Aj,  x1  and  a2.  and  to  £(g),  for  &  *  10  ,  10  ,  10  ,  1,  10, 

2  5 

10  ,  10  .  We  observe  that  the  distributions  of  all  the  systems  in  a  random 

environment  are  close  to  the  distribution  for  the  system  with  a  unique, 

average  rate  x^,  with  increasing  differences  when  the  environment  changes 

-5 

more  slowly.  We  also  observe  that  C(10  )  *  Y^X^)  +  Y2H(x2)»  uP  t0  five 

decimal  places. 

We  conclude  therefore  that  for  the  present  model,  the  random  environment 
has  little  effect  on  the  marginal  stationary  distribution  of  the  number  of 
machines  on  repair. 

In  order  to  measure  the  influence  of  the  random  environment  on  the 
dynamic  behaviour  of  the  system,  we  have  computed  the  average  time  needed  to 
-each  the  states  [n  machines  on  repair],  1  <  n  <  5,  starting  from  the  state 
[d  machine  on  repair].  We  present  on  Figure  1  the  results  for  n=5,  which  is 
terpretable  as  the  "time  to  complete  failure".  The  functions  are  as  follows 
1  and  2  :  f^(3)  and  f g ( 3 ) >  where  fj(0)  3  Ettime  to  reach  the  state  tall 
machines  on  repair]  starting  from  the  state  [0  machine  on  repair  and 
environment  j] }. 

-  ?,  4  and  5  :  g^),  g( x2 )  and  g(x^),  where  g(X)  *  E{time  to  reach  the  state 
[all  machines  on  repair],  starting  from  [0  machine  on  repair]}, for  a  system 
with  constant  failure  rate 
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-  6  :  y^(S)  +  Y2^2^)»  ecJua^  t0  the  stationary  expected  time  to  complete 
failure,  starting  from  [0  machine  on  repair]. 

We  clearly  observe  several  typical  ranges  of  values  for  B.  If  b  >  10, 
the  environment  changes  so  rapidly  that  the  expected  time  to  complete  failure 
is  equal  to  the  time  for  a  system  with  unique  failure  rate  equal  to  xQ. 

For  10'1<  6<10,  the  expected  time  to  complete  failure  does  not  depend  on  the 

r- 

initial  environment  state.  For  6  <  10“°,  the  environment  changes  so  slowly 
that  complete  failure  is  reached  before  a  change  of  environment  occurs. 

In  fact,  it  seems  that  the  system  is  almost  completely  decomposed  in  two 
different  systems,  one  corresponding  to  each  environment,  with  very  slow 

_3 

migrations  from  one  to  the  other.  Also,  we  note  that  for  &  <  10  ,  the  sta¬ 

tionary  expected  time  to  complete  failure  looses  any  practical  significance. 

Finally,  the  average  failure  rate  xQ  yields  an  overestimation  of  the  time  to 

-3 

complete  failure  for  10  <&<1,  that  is,  for  values  of  B  which  are  neither  much 
lower,  nor  much  higher  than  the  failure  rates. 

6.2  Mass  search  model .  We  define  a  reference  model  with  p=15  prey  and  m 
predators.  The  system  can  be  in  K=2  environments,  with  parameters  and 
c*2  =  2  aj.  The  detection  and  loss  rates  are  respectively  given  by  x^  =  .001, 
x^  =  .005,  m  =  .02  and  =  .01.  Thus,  environment  2  is  more  favorable  to  the 
predators,  since  the  detection  rate  is  higher,  and  the  loss  rate  is  lower; 
however,  environment  1  lasts  on  the  average  longer  than  environment  2. 

For  this  model,  the  rate  of  changes  in  the  environment  has  an  influence 
on  the  marginal  stationary  distribution  of  the  number  of  prey  under  surveillance. 
In  Table  II,  we  give  that  distribution  for  five  different  models.  Columns  1 
and  2  correspond  to  models  with  a  unique  environment  (respectively  environments 
1  and  2) . 
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■j i Linir'.s  3  and  4  correspond  to  the  reference  model,  respectively  with 

-  Id*4  (slow  environment  changes),  and  =  1  (fast  environment  changes). 

ciumn  5  corresponds  to  a  model  with  a  unique,  average  environment  : 

2  1  2  1 

4  A ^  +  -y-  ^2*  u  =  ~y  +  —y  W2*  ,Niote  that  the  distribution  in  column 

has  two  modes,  respectively  corresponding  to  the  modes  in  columns  1  and  2. 

We  have  also  measured  how  the  system  responds  to  increases  in  the  effi- 

-y  of  the  predators.  We  have  considered  two  ways  for  the  predators  to 

ecome  more  efficient.  The  first  one  is  by  becoming  more  numerous,  the  second 

increasing  the  probability  of  detecting  a  free  prey,  or  by  decreasing  the 

^ability  of  losing  a  prey  under  surveillance.  We  present  on  Figures  2  and 

'  e  values  of  t  .,  n=5,lQ,15,  j=l,2,  where 

n.J 

t  •  =  E  [time  until  n  prey  are  under  surveillance  (  at  time  0, 
n ,  j 

0  prey  under  surveillance,  environment  state  is  j], 

tue  reference  model  with  =  .0001,  and  m  equals  15  to  45.  It  clearly 
i.;,  ,;rs  that  environment  2  is  more  favorable  for  the  predators.  The  curve  for 
1:  .  presents  a  plateau  which  we  explain  as  follows.  We  denote  by  T*  .  the 

'  j  ■  ** 

'me  until  n  prey  are  under  surveillance  given  that  initially  zero  prey  are 
'dor  surveillance,  for  a  system  with  a  unique  environment,  identical  to 
. ^I'onment  j.  For  values  of  m  less  than  32,  say,  T15,l  is  so  large  (greater 

C 

*-  iO  on  the  average)  that  the  reference  model  with  =  .0001  switches  to 

.rorment  2  before  15  prey  are  under  surveillance  (the  switch  occurs  on  the 

verage  after  104  units  of  time).  Once  the  model  is  in  environment  2,  it  takes 

:  f  iff  average  10  to  10  units  of  time  to  have  all  preys  under  surveillance. 

nis  occurs  before  the  environment  switches  back  to  1,  and  the  total  elapsed 

4 

e  is,  on  the  average,  approximately  10  . 
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On  Figures  4  to  7,  we  present  the  values  of  t  .,  n=5,10,15,  j*l,2, 

n » j 

for  models  derived  from  the  reference  model.  We  consider  =  .0001, 
m  equals  to  15  and  20,  and  we  multiply  the  rates  A.  by  v*T  ,  we  divide  the 

J 

rates  w  ■  by  vT"  •  with  y  >  1.  so  that  the  probability  ratio's  "probability  of 

J 

detecting  /  probability  of  losing"  are  uniformly  multiplied  by  y  in  each  state. 
These  figures  may  be  used  in  conjunction  with  Figures  2  and  3  to  measure  trade¬ 
off  such  as  the  following  one.  Suppose  we  start  from  the  reference  model  with 
m=15,  and  that  we  double  the  number  of  predators  to  m=30.  This  will  entail  a 
reduction  on  t15  which  can  be  measured  on  Figure  2.  We  can  then  measure  on 
Figure  4  that  the  probability  ratio  must  be  multiplied  by  a  factor  y  =3  in  order 
to  obtain  the  same  reduction. 

In  Table  III,  we  give  the  approximate  values  of  y  which  give  the  same  reduction 
as  doubling  the  number  of  predators,  for  m=15  and  m=20. 
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APPENDIX  A. 


Algorithm  B  :  evaluation  of  the  stationary  probability  distribution,  the 
first  and  second  moments  of  the  time  to  reach  level  m,  1  <  m  <  N,  starting 
from  level  0. 

B1  First  passage  times  from  level  0  to  higher  levels. 

Bl.l  Initial  values. 


(  (-C"1)  -  -(A*0))-1  ; 

(Equation  (2.2)) 

u{1)  (-C-1)  e;  u^u^  ; 

(3.12) 

v^  .  2(-C'1)  u{1);  v*1*  ; 

(3.18) 

G(1)(0)  -  (-Cq1)  A<0);  G(1>1)(0)  -  G(1)(0); 

(3.11) 

Ud)  .  (-O'1)2  A<°);  U^1*1)  <-  } 

(3.15) 

B1.2  Levels  2  to  N. 

for  n  =  2  to  N  do 

{  (-c'^)  «-  -(A(n_1)+  M(n_1^ (-c"}2)  A(n"2h_1  ; 

(2.3) 

u<">.  (-c;!l)(e  +  ^("-DuC"-1)  ; 

(3.13) 

u(^»nJ  +.  u(l*n_l)+  G(l.n-1)(0)  u(n)  ; 

(3.14) 

v(n)  <-  2(-C“li)(I  +  M^^W0'1))  u<"> 

+  (-C“^)  ^n_1Mn‘1); 

(3.19) 

(3.20) 


d.n) 

+-  v(l»n*l)+  2  G^,n”^-(0)  v^; 

(3.20) 

(n)(0) 

-  (-C'.1!)  A^-1*; 

(3.11) 

(1,n)(o) 

<-  G(1,n‘1)(0)  G^n^(0); 

(3.8) 

(">  «. 

(-C^j)(I  +  M(n_1)U(n'1))(-C'l1)  A(n_1); 

(3.16) 

d.n)  , 

L)(l,n-l)G(n)(0^l,n-l)(0)  y(n)  } 

(3.17) 

B2  Stationary  distribution. 

B2.1  Initialization  of  the  recurrence. 

{  CN  «■  A(N)+  A(N_1);  (2.3) 

TT^|  •<-  solution  of  =  0  ,  e  =  1  }  (2.5) 


32.2  Determination  of  ,  0  <  n  <  N. 


for  n  =  N-l  to  0  do 


{  In 


Ifi+1 


M(n+i) (_c-i) . 


(re-normal ize) 
N 

a  I  e  ; 

k=n  “*  " 


for  k  =  n  to  N  do 


t  ^  <-  o’^>  > 


(2.6) 


WPUPP— JP)— M  WHWI-W  ■  "III  r<r-. • .  —  * 

1 

i 

| 

i 

j 

TABLES 


Table  I.  Cumulative  probabilities-number  of  machines  on  repair. 


Table  II.  Marginal  distribution-number  of  prey  under  surveillance. 
Missing  numbers  are  less  than  5.10"®. 


Table  III.  Values  of  y  giving  the  same  reduction  as  doubling  the 
number  of  predators. 


0 

1 

2 

3 

4 

. 5639b 

.84593 

.95872 

.99256 

.99932 

.56395 

.84593 

.95872 

.99256 

.99932 

.56398 

.84589 

.95870 

.99255 

.99932 

.56421 

.84561 

.95850 

.99249 

.99931 

.56579 

.84397 

.95711 

.99202 

.99925 

.56910 

.84151 

.95424 

.99091 

.99908 

.57033 

.84078 

.95321 

.99047 

.99900 

.57050 

.84068 

.95306 

.99040 

.99899 

.49516 

.79226 

.93486 

.98620 

.99852 

.72118 

.93754 

.98946 

.99881 

.99993 

Table  I 


Cumulative  probabilities-number  of  machines  on  repair 


n 

I 

II 

III 

IV 

V 

0 

.00073 

.00048 

1 

.00820 

.00543 

.00005 

.00003 

2 

.04019 

.02662 

.00061 

.00043 

3 

.11319 

.07504 

.00428 

.00342 

4 

1 

.20375 

.00002 

.13524 

.01959 

.01724 

5 

.24653 

.00026 

. 16400 

.06175 

.05835 

6 

.20544 

.00220 

1 

.13776 

.13712 

.13601 

7 

.11886 

.01270 

.08407 

.21606 

.22012 

8 

.04755 

.05082 

.04948 

.24021 

.24628 

9 

.01294 

.13835 

.05546 

. 18523 

.18752 

10 

.00233 

.24903 

.08481 

.09610 

.09441 

11 

.00026 

.28299 

.09427 

.03198 

.03001 

12 

.00002 

.18866 

.06255 

.00633 

.00560 

13 

.06531 

.02160 

.00066 

.00054 

14 

15 

.00933 

.00031 

.00308 

.00010 

.00003 

.00002 

Table  II 

Limiting  distribution-number  of  prey  under  surveillance. 
Missing  numbers  are  less  than  5.10"® 


i 


I 


m 

*5,1 

t10,l 

t15 ,1 

t5 ,2 

_ 

t10,2 

t15,2 

15 

4.1 

4.4 

7.8 

D 

7.1 

B 

20 

4.1 

4.0 

3.9 

B 

5.6 

1 

Table  III 


Values  of  y  giving  the  same  reduction  as  doubling 
the  number  of  predators 


FIGURES 


Expected  time  to  reach  the  state  [all  machines  on  repair], 
starting  from  [no  machine  on  repair]. 

Expected  time  to  reach  the  states  [n  prey  under  observation], 
starting  from  [0  prey  under  observation,  environment  1], 
varying  number  of  predators. 

Expected  time  to  reach  the  states  [n  prey  under  observation], 
starting  from  [0  prey  under  observation,  environment  2], 
varying  number  of  predators. 

Expected  time  to  reach  the  states  [n  prey  under  observation], 
starting  from  [0  prey  under  observation,  environment  1], 

15  predators,  varying  probability  ratio. 

Expected  time  to  reach  the  states  [n  prey  under  observation], 
starting  from  [0  prey  under  observation,  environment  2], 

15  predators,  varying  probability  ratio. 

Expected  time  to  reach  the  states  [n  prey  under  observation], 
starting  from  [0  prey  under  observation,  environment  1], 

20  predators,  varying  probability  ratio. 

Expected  time  to  reach  the  states  [n  prey  under  observation], 
starting  from  [0  prey  under  observation,  environment  2], 

20  predators,  varying  probability  ratio. 


Fig. 4.  Expected  time  to  reach  the  states  [n  prey  i 
from  [0  prey  under  observation,  environmen 
probability  ratio. 
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